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Improving the triaxial bulge model of M31 



Simon Berman 1 and Laurent Loinard 2 



(N 
O 
O 
(N 

C 



> 

in 
oo 
(N 

O 
(N 
O 

Of 

6 



1 Theoretical Physics, University of Oxford, 1 Keble Road, Oxford, UK; simon@thphys.ox.ac.uk 

2 Instituto de Astronomia, Universidad Nacional Autonoma de Mexico, Apartado Postal 72-3, 58089 Morelia, Michoacdn, Mexico 



Accepted . Received ; in original form 



ABSTRACT 

A detailed hydro-dynamical model of the gas flow in the triaxial gravitational potential 
of the bulge of the Andromeda galaxy (M31) has recently been proposed by Berman 
(2001), and shown to provide excellent agreement with the CO emission line velocities 
observed along its major axis. In the present paper, we confirm the validity of that 
model by showing that it can also reproduce the CO velocities observed off the major 
axis - a much more robust test. The CO observations, however, tend to span a wider 
range of velocities than a direct application of the original model of Berman would 
suggest. This situation can be improved significantly if the molecular disk is made 
thicker, a requirement already encountered in dynamical simulations of other spiral 
galaxies, and typically attributed to a broadening of the molecular layer in galactic 
fountain-like processes. In the central regions of M31, however, it is unclear whether 
there actually is a thick molecular disk, or whether broadening the molecular layer 
is merely an artificial theoretical means of accounting for some disk warping. Other 
effects not included in the model, such as hydraulic jumps, might also contribute to a 
widening of the velocities. 

Key words: hydrodynamics - galaxies: M31 - galaxies: ISM - galaxies: structure - 
galaxies: evolution. 



1 INTRODUCTION 

There is growing evidence that spiral galaxies in which no 
bar is easily seen on optical or infrared images, nonethe- 
less, have triaxial bulges. Twisting of the inner isophotes 
and misalignments between the disk and the bulge major 
axes (both of which have been detected in several spirals 
with no obvious bar) are clear signposts of triaxiality. Even 
stronger cases have been made by combining optical or in- 
frared photometric data with spectroscopic observations of 
the gaseous interstellar component. Indeed, when submit- 
ted to a triaxial gravitational potential, the interstellar gas 
near the centre of a spiral galaxy can be found at higher 
velocities than expected from circular motion. Such anoma- 
lous velocities have been reported in several galaxies where 
no strong bar is visible, including the two nearest examples: 
the Milky Way (Rougoor & Oort 1959; Dame et al. 2001) 
and the Andromeda galaxy (Lindblad 1956; Loinard et al. 
1999). 

The triaxiality of spiral galaxy bulges is particularly im- 
portant for two reasons. First, it can affect the dynamics of 
interstellar gas, gathering large amounts of it near galactic 
centres. This could help to explain the appearance of active 
galactic nuclei. On long timescales, a build of gas in galactic 
centres could affect the morphology of spiral galaxies, and 
imply secular changes in their Hubble types. Second, triax- 



iality allows more constraints to be derived from spectral 
observations. In an axisymmetric system, the rotation curve 
only allows the determination of the radial mass distribu- 
tion, without providing any three dimensional information. 
In particular, it is not possible to assess in an axisymmetric 
system whether the dark matter component is distributed 
in a flat disk-like structure or in a spherical halo. In a triax- 
ial system, that disk-halo degeneracy can be removed, with 
obviously important consequences for our understanding of 
dark matter. Cold Dark Matter (CDM) models, for instance, 
favor spherical halos. 

Since triaxiality primarily affects the dynamics at the 
centres of spiral galaxies, it should be sought in our neigh- 
bours, where the achievable linear resolution is highest. This 
is particularly important when using radio observations of 
the gaseous component, which often have limited angular 
resolutions. Detailed studies of the dynamics of the cen- 
tral regions of the Milky Way are adversely affected by our 
location inside of the system, and the well-known associ- 
ated confusion and distance ambiguity problems. However, 
the existence of triaxial features in the inner Galaxy have 
been postulated for over forty years. Based on the discovery 
of the 3 kpc arm, Kerr (1967) proposed the existence of a 
gaseous bar at the Galaxy centre, and Mulder & Liem (1986) 
modelled the arm as a density wave in a barred potential. 
By combining photometric and kinematic data, Gerhard & 
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Figure 1. Comparison of CO observations and gas velocities from the model of maximum likelihood of Berman (2001) away from the 
line of nodes of the disk. Model data were convolved with the 2D function representing the FCRAO beam as described in the text. 
Intensity is plotted as a greyscale for observations and using contour levels for model data. Contours are given at 1, 3, 10 and 30% of 
the maximum value in each subplot. The dashed line indicates the systemic velocity and the solid line the results at Y arcsec ignoring 
the effects of the thick disk and of the beam function. 



Vietri (1986) demonstrated the triaxiality of the Galactic 
bulge, and COBE images of the Galactic centre have con- 
firmed this triaxiality (Dwek et al. 1995; Binney et al. 1997). 

The next nearest spiral (the Andromeda Galaxy - M31) 
has no clear visible bar, and is an excellent target for a dy- 
namics study, because its high inclination enables an accu- 
rate determination of the gas kinematics. M31 has long been 
known to have twisted inner isophotes and misaligned bulge 
and disk major axes (Lindblad 1956). Moreover, the exis- 
tence of anomalous gas velocities in the inner few kilopar- 
secs is also now well established both for the atomic (Brinks 
& Burton 1984) and molecular (Loinard et al. 1996, 1999) 
components. A convincing interpretation in terms of triax- 
iality was put forward as early as 1956 by Lindblad, and 
was further developed by Stark (1977) and Stark & Binney 
(1994). However, those early works were hampered by the 
lack of strong observational constraints. 

The distribution and kinematics of the atomic compo- 
nent of the interstellar medium (ISM) has long been known 
at high resolution and sensitivity across the entire disk of 
M31 thanks to interferometric observations of the 21-cm line 
of atomic hydrogen (Brinks & Shane 1984, and references 
therein). However, those observations do not provide clear 
information on the ISM kinematics in central regions, be- 
cause the kinematic component associated with the warped 
outer disk can be seen in projection through the inner disk. 
The resulting confusion between inner and outer disk seen 



in the H I position- velocity diagrams of the central regions of 
M31 precludes accurate studies of its dynamics. The molec- 
ular component traced by CO emission is a better tracer of 
the dynamics of the inner regions because no CO emission 
can be detected in the warped outer disk. However, the in- 
ner regions of M31 are also particularly dim in CO (Dame 
et al. 1993, Loinard et al. 1999). Fragmentary CO observa- 
tions of the central regions of M31 were used to constrain the 
model of Stark & Binney (1994). Recently, however, more 
systematic CO observations have been obtained. 

Loinard et al. (1995) presented a deep search for CO 
emission along the inner major axis of M31. This search con- 
firmed the existence of anomalous velocities in the inner few 
kiloparsecs of M31, and results along the line of nodes of the 
disk were used by Berman (2001) to constrain an improved 
dynamical model of the central regions of M31 (see §2 be- 
low). However, Berman (2001) pointed out that off-axis CO 
observations would provide better constraints. Such data are 
now available - at least for the southern half of M31 - thanks 
to the CO(1-0) survey made at the Five College Radio As- 
tronomy Observatory (FCRAO - Loinard et al. 1996, 1999). 
The angular resolution of the survey is ~ 1 arcmin, and it 
covers the entire Southern part of M31 with a sampling of 
50 arcsec. The observational noise level is rather constant 
across the entire surveyed region with a typical r.m.s. of 
45-50 mK per 3.25 km s _1 spectral channel. Although the 
CO (1-0) is found to peak in the broad Population I ring at 
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Figure 2. Comparison of CO observations and gas velocities from the model of maximum likelihood away from the line of nodes of 
the disk. The model includes a corotating thick disk of scale height zq = 400 pc. Model data were convolved with the 2D function 
representing the FCRAO beam as described in the text. Intensity is plotted as a greyscale for observations and using contour levels 
for model data. Contours are given at 1, 3, 10 and 30% of the maximum value in each subplot. The dashed line indicates the systemic 
velocity and the solid line the results at Y arcsec ignoring the effects of the thick disk and of the beam function. 



10 kpc from the centre, significant emission is still detected 
in the inner 2-3 kpc (10-15 arcmin), where most constraints 
can be obtained about the triaxiality of the bulge. The inte- 
grated intensity images or position-velocity diagrams shown 
in Loinard et al. (1999), as well as the position-velocity dia- 
grams that will be shown in this article, have been obtained 
after 'unsharp masking' was applied to the data. In this 
scheme, a data cube smoothed both in position and veloc- 
ity (but not regridded) is first constructed from the origi- 
nal data. The average noise level a in this smoothed cube 
is computed, and all the pixels in the non-smoothed data 
cube where no significant (> 3a) emission is found in the 
smoothed version are blanked. This avoids that unnecessary 
noise is added when summations over several pixels are per- 
formed, and in effect reduces the noise in integrated maps 
by a factor of a few. A slight drawback of this method is 
that the noise becomes dependent on the number of pixels 
summed (or averaged), and, therefore, on the position in the 
maps. 



Anomalous velocities are clearly seen in this data set 
in the inner regions of M31 (see fig. 12 in Loinard et al. 
1999). While the cut along the major axis essentially shows 
the same anomalous velocities as reported by Loinard et 
al. (1995), cuts parallel to the major axis provide new con- 
straints for the models, as called for by Berman (2001). 



2 ORIGINAL MODEL 

The dynamical model developed by Berman (2001) takes 
into account a rotating, triaxial bulge, and an axisymmetric 
component that mimics the combined influence of the disk 
and halo. The properties of the bulge are obtained by finding 
a triaxial mass distribution which gives the observed R 1 ^ 4 
B-band surface brightness profile when rotated and inclined 
to the plane of the sky. Properties of the axisymmetric com- 
ponent are deduced from a simple interpretation of H I and 
CO kinematics. The gas response to the gravitational poten- 
tial is computed using a hydrodynamical code (galahad) 
that solves the non-self gravitating, isothermal Euler equa- 
tions, on a regular Cartesian grid, using the FS2 algorithm 
of van Albada et al. (1982). Each grid cell is 125 pc (37.5 
arcsec) on a side, which is similar to the resolution (~ 60 
arcsec) of the CO data. Removal and injection of interstellar 
matter by star-formation and stellar mass-loss respectively 
are included (somewhat crudely) in the algorithm. 

Berman (2001) constrained the parameters of his model 
by comparing the output gas flow to the CO observations 
along the line of nodes of the disk of M31 of Loinard et al. 
(1995). The best fit to the data is obtained for a bulge semi- 
major axis a = 3.5 kpc, a pattern speed fl p — 53.7 km s _1 
kpc -1 , a corotation radius to bulge semi-major axis ratio 
TZ = 1.2, an angle <f> between the major axis of the bulge 
and the disk's line of nodes of 15°, and a bulge mass-to- 
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Figure 3. Comparison of CO observations and gas velocities from the model of maximum likelihood away from the line of nodes of the 
disk. The model includes a thick disk of scale height zq = 200 pc and a loss of velocity of 0.3 vo kpc -1 away from the plane of the 
disk. Model data were convolved with the 2D function representing the FCRAO beam as described in the text. Intensity is plotted as 
a greyscale for observations and using contour levels for model data. Contours are given at 1, 3, 10 and 30% of the maximum value in 
each subplot. The dashed line indicates the systemic velocity and the solid line the results at Y arcsec ignoring the effects of the thick 
disk and of the beam function. 



light ratio Tb = 6.5. This in turns implies a bulge massy 
of M bulgc = 2.3 x 10 10 Mq . The disk mass within the ra- 
dius of the bulge (3.5 kpc) can be estimated by assuming 
reasonable values for the B-band exponential disk central 
surface brightness (Jo = 21.6 mag arcsec -1 , Walterbos & 
Kennicutt 1988) which when corrected for absorption and 
inclination using the values of Berman (2001) gives a cen- 
tral surface density of Ed =113 L© pc -3 , disk scale length 
(Rd = 5.8 kpc, Walterbos & Kennicutt 1988), and a conser- 
vative disk mass-to- light ratio of 4. This gives a disk mass 
of Mdisk = T2 x 10 10 M . The total (bulge + disk + halo) 
mass inside the same radius can be estimated from the ro- 
tation curve as M tota i = 3.7 xlO 10 Mq, leaving little room 
for a dark halo. 

The value of 1Z — 1.2 found by Berman (2001) implies 
that the bulge of M31 is a fast rotator. Interestingly, of the 
few galaxies for which 1Z has been measured, all have been 
found to be fast rotators. Debattista & Sellwood (2000) 
showed that such fast bars cannot co-exist with massive 
dark halos because dynamical friction would slow them down 
rapidly. The value of 1Z = 1.2 is further evidence that M31 
must have a minimal halo, and, consequently, a maximum 



1 The value given for the mass of the bulge in Berman (2001) was 
not correct. The amended version has been used in this paper. 



disk. This is notably at odds with the contentions of Cold 
Dark Matter (CDM) models, which predict massive halos 
and minimal disks (e.g. Navarro et al. 1997). 



3 MODELLING OFF-AXIS OBSERVATIONS 

As mentioned earlier, Berman (2001) pointed out that a 
much stronger case could be made for his best fit model 
if it were found to be consistent with observations of CO 
velocities away from the line of nodes of the disk. Moreover, 
it is important to utilize these off-axis velocities since dy- 
namical models of the bulge of M31 have never before been 
constrained by real 3D (X, Y, V) observations. 

To compare accurately the output of the model with 
the data, we performed a 2D convolution of the model with 
a realistic description of the FCRAO beam, consisting of 
the sum of a sine function and two Gaussians. The side- 
lobes of the sine function provide a good (very conservative) 
description of any coma effects, while the two Gaussians 
describe the response of the FCRAO antenna to emission 
coming from angles far away from the nominal pointing po- 
sition, due to the imperfections of the dish (Ladd &i Heyer 
1996). The sine function has a full width at half maximum 
(FWHM) of 45 arcsec and the two Gaussians have FWHM 
and attenuations of 30 arcmin and 30 dB, and 4 degrees and 
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50 dB, respectively. This convolution produces multiple ve- 
locities at each X position in the (X, V) position-velocity 
diagrams derived from the model, in the same way as mul- 
tiple velocities arise naturally in the data cube. 

Fig.^ shows position-velocity diagrams for observa- 
tional data at Y — ±y x 50 arcsec, where y = 1, 2, 3. Over- 
layed are Berman's best fit model's predictions convolved 
with the FCRAO beam function which should, if the model 
is correct, approximately enclose the observations on every 
plot. In this, as in all of the comparisons below, the model 
has been interpolated to fit the positions of the observa- 
tions. Although the general trends of the predictions and 
the observations are similar, the model data enclose a nar- 
rower range of velocities than the observations and are fur- 
ther from the systemic velocity than the observations. The 
discrepancies are just as large or greater for models with dif- 
ferent input parameters. In particular, it is not reasonable to 
use a smaller rotation velocity. This would help reconcile the 
predictions of the model proposed by Berman (2001) with 
the data in the central regions, but would predict velocities 
much too small at distances larger than about 40 arcmin 
from the centre. 

The same 'beard' phenomenon is seen in the position- 
velocity diagrams of the spiral galaxy NGC 2403 (fig. 1 of 
Schaap et al. 2000). The arguments presented in that paper 
conclude that the beard is either caused by a thick high den- 
sity gas layer with a half width at half maximum (HWHM) 
intensity of 500 pc or a lower density layer with HWHM of 
1.75 kpc, rotating more slowly than the standard thin disk. 
This conclusion is motivated by galactic fountain models 
(Spitzer 1990) in which hot gas from supernova explosions 
and galactic winds rises from the disk into the halo. As it 
rises, the gravitational attraction towards the centre of the 
galaxy lessens, the gas moves outwards and, by conserva- 
tion of angular momentum, its azimuthal velocity decreases. 
Thick molecular disks have been seen in CO in the Milky 
Way Galaxy (Dame & Thaddeus 1994), and the edge-on 
spiral galaxy NGC 891 (Garcia-Burillo et al. 1992). 

To recover the small magnitude and large range of veloc- 
ities of the observations of M31, the bulge model of Berman 
(2001) is augmented by the kinematics of a thick disk which 
has been rotated and inclined to the observer's frame. The 
vertical density distribution follows Schaap et al. (2000) in 
taking a Gaussian form 

p(z) = po exp(-z 2 /2z%), (1) 

P ° = (27^/2 V 

where po and So are the density and surface density respec- 
tively at z — 0, and zo is the vertical scale height of the gas. 
The half width at half maximum intensity Z]/ 2 = 1-18 Zo- 
The intensities at (X, Y) are convolved in 2D with the func- 
tion representing the FCRAO beam as described earlier, to 
produce position-velocity plots for Y = ± y x 50 arcsec, 
where y = 1,2,3. Fig. | compares the observational data 
with the results from a model with zo = 400 pc, the lowest 
value of the scale height that adequately encloses the vast 
majority of the observations for X > — 50 arcmin. For this 
model, Zi/2 = 472 pc. This is very large when compared 
to the thick disk observed by Dame & Thaddeus (1994) in 
the Milky Way, for which HWHM is between 71 and 133 
pc depending on position, but smaller than the thick disk 



of NGC 891 (Garcia-Burillo et al. 1992), where HWHM is 
larger than 1 kpc. 

We now follow Schaap (2000) and drop the implicit as- 
sumption that the thick disk is corotating with the z = 
plane and reduce the velocities of the gas where z 7^ 0. Ve- 
locities parallel to the z axis stay fixed at zero. A slow thick 
disk is implied by the galactic fountain models and was also 
invoked by Swaters et al. (1997) to explain the Hi observa- 
tions of NGC 891. We take the linear velocity distribution 
used in both papers: 

v(z) = v (l-7«), (3) 

where vo is the velocity in the z = plane, 7 is the velocity 
lost per kpc and z is the height. By slowing the gas down 
above and below the plane of the disk, the range of velocities 
at a particular scale height is extended. It can be seen in Fig. 
|§] that a similar range of velocities to both the observations 
and the corotating model is produced if the scale height 
is just 200 pc and 7 = 0.3 i.e. the gas loses 30% of its 
speed per kpc away from the plane of the disk. This implies 
that Z1/2 = 236 pc which is closer to but still larger than 
the values quoted for the Milky Way in Dame & Thaddeus 
(1994) and much smaller than the value for NGC 891. 

By including a thick disk, the model initially proposed 
by Berman (2001) is able to match the 3D CO data of 
Loinard et al. (1999) reasonably well. The most notice- 
able discrepancy is found around X — ±20 — 25 arcmin, 
V ~ —45 km s _1 , a position that the original model found 
difficult to account for even along the line of nodes of the 
disk. It is plausible that some vertical or turbulent gas mo- 
tions not included in our model can account for these re- 
maining discrepancies. 

As mentioned above, galactic fountain models are often 
invoked to account for thick gaseous disks in the centres of 
galaxies. However, they are probably not justifiable in M31 
since it does not harbour a great deal of star formation. 
Moreover, whilst the thick disk found by Dame & Thad- 
deus (1994) in the Milky Way appears to be much fainter 
in CO than the Galactic mid-plane, the present modelling 
requires that in M31 it be just as bright. A plausible al- 
ternative follows from the recent work by Martos & Cox 
(1998), who describe many numerical simulations combin- 
ing magnetic fields with hydrodynamics to create postshock 
regions in which large pressure gradients exist near to the 
mid-plane causing the gas layer to expand vertically. The 
result is a 'hydraulic jump' in which a build up of gas in the 
vertical direction can be seen with scale heights of hundreds 
of parsecs. 

Another, perhaps more controversial, possibility should 
be mentioned here. Pringle et al. (2001) contend that much 
of the ISM in spiral galaxies is molecular, but too cold to be 
easily detectable in CO. In this view, only the parts of molec- 
ular clouds illuminated by newly formed stars are bright in 
CO. The CO-bright thin disk could then correspond to the 
thin layer where stars form most actively, whereas the thick 
disk would correspond to the real total extent of the molec- 
ular layer. The existence of large quantities of cold molec- 
ular gas in the inner regions of M31 has been proposed by 
Loinard & Allen (1998 - and references therein). Moreover, 
since the inner region of M31 harbours little star-forming ac- 
tivity, the thin disk there would not be particularly brighter 
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in CO than the thick disk, and the distinction between the 
two would indeed vanish. 

Finally, a last possibility is that the inner disk of M31 is 
actually thin, but warped. Modelling it as a thick disk then 
becomes a convenient way to account for the gas located in 
the warped parts of the disk, above or below the mid-plane. 



4 CONCLUSIONS 

In this article, we show that hydrodynamical triaxial mod- 
els of the bulge of the Andromeda galaxy (M31) previously 
presented by Berman (2001), and initially tested against CO 
emission observations along the apparent major axis only, 
can also account for the CO emission velocities observed off 
the apparent major axis - a much more robust validity test. 
The triaxial model implies that the bulge of M31 is a fast 
rotator and hence the dark matter contained in its central 
regions must be minimal and coincident with the stellar disk. 

To account for the whole velocity range covered by the 
observations, a finite, fairly large thickness (half width at 
half maximum intensity of 200-500 pc) must be given to 
the molecular disk. It is quite unclear, however, whether the 
inner molecular disk of M31 is truly that thick or whether 
adding it to the model is merely an effective theoretical trick 
able to account for a warping of a thiner disk. Other out of 
plane, or peculiar velocities (such as those associated with 
hydraulic jumps or turbulence) could also contribute to the 
wide velocity range covered by the CO observations. 
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